A synthesis of dryland restoration techniques.

Purpose

To quantitatively examine the efficacy of vegetation restoration in drylands globally.

Questions

  1. What is the global extent of research that directly examined restoration of drylands?
  2. What were the common measures?
  3. Is the restoration of vegetation a common and primary focus?
  4. How frequently does the restoration measure outcomes beyond the focal species?
  5. What were the primary restoration goals as reported by primary authors?
  6. How much variation was there in the techniques tested and how long were experiments monitored and tested?
  7. How relatively effective were the techniques?

Step 2. Sort

A summary of sort process using PRISMA.

library(PRISMAstatement)
prisma(found = 1504,
       found_other = 5,
       no_dupes = 1039, 
       screened = 1039, 
       screen_exclusions = 861, 
       full_text = 178,
       full_text_exclusions = 101, 
       qualitative = 77, 
       quantitative = 66,
       width = 800, height = 800)

Step 3. Synthesize

Check data and calculate necessary measures.

#all data
#data <- read_csv("data/data.csv")
#data <- data %>%
  #mutate(lrr = log(mean.t/mean.c), rii = ((mean.t-mean.c)/(mean.t + mean.c)), var.es = ((sd.t^2/n.t*mean.t^2) + (sd.c^2/n.c*mean.c^2)))

#other effect size estimates
#library(compute.es)
#data <- data %>%
  #mutate(d=mes(mean.t, mean.c, sd.t, sd.c, n.t, n.c, level = 95, #cer = 0.2, dig = 2, , id = ID, data = data))
#check metafor

#data from ag and grazing studies that examined restoration in drylands
data <- read_csv("data/data.csv")
mydata<- data %>% 
 filter(disturbance == c("agriculture","grazing")) %>% 
 filter(!notes %in% "couldnt extract data")


mydata <- mydata %>%
  mutate(lrr = log(mean.t/mean.c), rii = ((mean.t-mean.c)/(mean.t + mean.c)), var.es = ((sd.t^2/n.t*mean.t^2) + (sd.c^2/n.c*mean.c^2))) #use only lrr now

Step 4. Summarize

Explore summary level data of all data. Explore aggregation levels that support the most reasonable data structure and minimize non-independence issues.

#evidence map####
require(maps)
world<-map_data("world")
map<-ggplot() + geom_polygon(data=world, fill="gray50", aes(x=long, y=lat, group=group))
#map + geom_point(data=paperdata, aes(x=long, y=lat)) +
  #labs(x = "longitude", y = "latitude") #render a literal map, i.e. evidence map, of where we study the niche in deserts globally

#add in levels and color code points on map####
#map of 178 articles
map + geom_point(data=data, aes(x=long, y=lat, color = paradigm)) + 
  scale_color_brewer(palette = "Paired") +
  labs(x = "longitude", y = "latitude", color = "")

#map of selected articles, agriculture and grazing disturbances
map + geom_point(data=mydata, aes(x=long, y=lat, color = paradigm)) + 
  scale_color_brewer(palette = "Paired") +
  labs(x = "longitude", y = "latitude", color = "")

#aggregation####
se <- function(x){
  sd(x)/sqrt(length(x))
}

data.simple <- mydata %>%
  group_by(study.ID, paradigm, technique, measure.success) %>%
  summarise(n = n(), mean.lrr = mean(lrr), mean.rii = mean(rii), mean.var = mean(var.es))

main.data <- mydata %>%
  group_by(study.ID, paradigm, intervention, outcome) %>%
  summarise(n = n(), mean.lrr = mean(lrr), mean.rii = mean(rii), mean.var = mean(var.es))


#EDA data####
simple.data <- mydata %>% group_by(study.ID, paradigm, technique, measure.success) %>% summarise(mean.rii = mean(rii), error = se(rii))
simple.data <- na.omit(simple.data)

parad.data <- mydata %>% group_by(study.ID, paradigm) %>% summarise(mean.rii = mean(rii), error = se(rii))
parad.data <- na.omit(parad.data)

tech.data <- mydata %>% group_by(study.ID, technique) %>% summarise(mean.rii = mean(rii), error = se(rii))
tech.data <- na.omit(tech.data)

success.data <- mydata %>% group_by(study.ID, measure.success) %>% summarise(mean.rii = mean(rii), error = se(rii))
success.data <- na.omit(success.data)


#active
active <- mydata %>%
  filter(paradigm == "active")

#viz for aggregation####
disturbance.data <- mydata %>% 
  group_by(study.ID,disturbance, paradigm) %>% 
  count()

disturbance.data2 <- disturbance.data %>% 
  group_by(disturbance,paradigm) %>% 
  count()

ggplot(na.omit(disturbance.data2), aes(disturbance,nn, fill=paradigm))+ 
  geom_bar(stat = "identity") + 
  coord_flip(ylim=0:44) + 
  scale_fill_brewer(palette = "Blues")

intervention.data <- active %>% 
  group_by(study.ID,intervention, paradigm) %>% 
  count()

intervention.data2 <- intervention.data %>% 
  group_by(intervention,paradigm) %>% 
  count()

#ggplot(na.omit(intervention.data2), aes(intervention,nn, fill=paradigm)) +
  #geom_bar(stat = "identity") + 
  #coord_flip(ylim=0:44) + 
  #scale_fill_brewer(palette = "Blues")

technique.data <- mydata %>% 
  group_by(study.ID,technique, paradigm) %>% 
  count()

technique.data2 <- technique.data %>% 
  group_by(technique,paradigm) %>% 
  count()

technique.data3 <- technique.data2 %>% 
  group_by(technique) %>% 
  count()


ggplot(na.omit(data.simple), aes(technique, n, fill = paradigm)) + 
  geom_bar(stat = "identity") + 
  coord_flip() + 
  scale_fill_brewer(palette = "Paired")

ggplot(na.omit(data.simple), aes(measure.success, n, fill = paradigm)) + 
  geom_bar(stat = "identity") + 
  coord_flip() + 
  scale_fill_brewer(palette = "Paired")

#better
ggplot(main.data, aes(intervention, n, fill = paradigm)) + 
  geom_bar(stat = "identity") + 
  coord_flip() + 
  scale_fill_brewer(palette = "Paired") +
  labs(fill = "")

ggplot(main.data, aes(outcome, n, fill = paradigm)) + 
  geom_bar(stat = "identity") + 
  coord_flip() + 
  scale_fill_brewer(palette = "Paired") +
  labs(fill = "")

Step 5a. EDA

Exploratory data analyses to understand data and QA/QC using Rii.

Step 5b. Models

Meta and conventional statistical models to explore relative efficacy.

Step 5. Synthesis stats

#p-value meta
#library(metap)
#all data for metas but cleaned for na's
#all data for meta cleaned
mdata.all <- mydata %>%
  filter(!is.na(lrr)) %>%
  filter(!is.na(var.es)) %>%
  filter(!is.na(n.t)) %>%
  filter(!is.na(p)) %>%
  filter(!is.na(intervention)) %>%
  filter(is.finite(lrr)) %>%
  rename(exp.length = `exp.length (months)`) #%>%
  #filter(!is.na(exp.length))

#active only data for meta cleaned up
mdata <- mydata %>%
  filter(paradigm == "active") %>%
  filter(!is.na(lrr)) %>%
  filter(!is.na(var.es)) %>%
  filter(!is.na(n.t)) %>%
  filter(!is.na(p)) %>%
  filter(!is.na(intervention)) %>%
  filter(is.finite(lrr)) %>%
  rename(exp.length = `exp.length (months)`) #%>%
  #filter(!is.na(exp.length))

#aggregated data for metas var estimated with central tendency
#note - could also bootstrap mean variance here instead of arithematic mean
simple.mdata <- mdata %>%
  group_by(intervention) %>%
  summarise(lrr = mean(lrr), var.es = mean(var.es), n = mean(n.t))

simple.mdata.2 <- mdata %>%
  group_by(intervention, outcome) %>%
  summarise(lrr = mean(lrr), var.es = mean(var.es), n = mean(n.t))  

simple.mdata.var <- mdata %>%
  group_by(intervention) %>%
  summarise(lrr = mean(lrr), var.es = se(var.es), n = mean(n.t))

simple.mdata2.var <- mdata %>%
  group_by(intervention, outcome) %>%
  summarise(lrr = mean(lrr), var.es = se(var.es), n = mean(n.t))

#metas with p-values####
#schweder(mdata$p)
#sumz(p, data = mdata)
#mdata %>%
  #split(.$intervention) %>%
  #purrr::map(~sumz(.$p)) 
#sumlog(mdata$p)

#metas with meta package on effect size measures####
library(meta)
#all data non-aggregated
#m <- metagen(lrr, var.es, studlab = ID, byvar = intervention, data = mdata)
#summary(m)
#funnel(m)
#metabias(m)
#forest(m, sortvar = intervention)

#t-tests if different from 0

tmu <- function(x){t.test(x, mu = 0, paired = FALSE, var.equal=FALSE, conf.level = 0.95)}

mdata.all %>%
  split(.$paradigm) %>%
  purrr::map(~tmu(.$lrr)) #note this uses arithmetic means not estimated means from random effect models
## $active
## 
##  One Sample t-test
## 
## data:  x
## t = 3.4074, df = 359, p-value = 0.0007302
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  0.1316402 0.4909925
## sample estimates:
## mean of x 
## 0.3113163 
## 
## 
## $passive
## 
##  One Sample t-test
## 
## data:  x
## t = -5.8407, df = 171, p-value = 2.559e-08
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.5034850 -0.2491308
## sample estimates:
##  mean of x 
## -0.3763079
#Propose we only use random effect models because of the diversity of studies
m <- metagen(lrr, var.es, studlab = ID, byvar = paradigm, comb.fixed=FALSE, data = mdata.all)
summary(m) #active is positive and passive is negative so should not mix
## Number of studies combined: k = 496
## 
##                                         95%-CI     z p-value
## Random effects model -0.0139 [-0.0855; 0.0578] -0.38  0.7045
## 
## Quantifying heterogeneity:
## tau^2 = 0.2648; H = 12688836.57 [12688835.55; 12688837.59]; I^2 = 100.0% [100.0%; 100.0%]
## 
## Test of heterogeneity:
##                     Q d.f. p-value
##  79698253852777280.00  495       0
## 
## Results for subgroups (random effects model):
##                      k                     95%-CI                    Q
## paradigm = active  333  0.1947 [ 0.1033;  0.2862] 79607575146295024.00
## paradigm = passive 163 -0.3537 [-0.4882; -0.2192]    90524299681441.44
##                     tau^2    I^2
## paradigm = active  0.2649 100.0%
## paradigm = passive 0.3666 100.0%
## 
## Test for subgroup differences (random effects model):
##                      Q d.f.  p-value
## Between groups   43.68    1 < 0.0001
## 
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
#metabias(m)
#should do a t-test right now of paradigm is diff from 0

#interventions
m <- metagen(lrr, var.es, studlab = ID, byvar = intervention, comb.fixed=FALSE, subset = paradigm == "active", data = mdata.all)
summary(m)
## Number of studies combined: k = 333
## 
##                                       95%-CI    z  p-value
## Random effects model 0.1947 [0.1033; 0.2862] 4.17 < 0.0001
## 
## Quantifying heterogeneity:
## tau^2 = 0.2649; H = 15484891.12 [15484889.86; 15484892.37]; I^2 = 100.0% [100.0%; 100.0%]
## 
## Test of heterogeneity:
##                     Q d.f. p-value
##  79607575146295024.00  332       0
## 
## Results for subgroups (random effects model):
##                                 k                   95%-CI
## intervention = vegetation     221 0.1403 [ 0.0378; 0.2428]
## intervention = soil            75 0.4617 [ 0.3633; 0.5601]
## intervention = water addition  37 0.0000 [-0.0014; 0.0014]
##                                                  Q  tau^2    I^2
## intervention = vegetation     79607575088050608.00 0.2649 100.0%
## intervention = soil                     2099026.81 0.0509 100.0%
## intervention = water addition                 2.34      0   0.0%
## 
## Test for subgroup differences (random effects model):
##                      Q d.f.  p-value
## Between groups   91.80    2 < 0.0001
## 
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
#metabias(m)

m <- metagen(lrr, var.es, studlab = ID, byvar = intervention, subset = paradigm == "passive", comb.fixed=FALSE, data = mdata.all)
summary(m)
## Number of studies combined: k = 163
## 
##                                          95%-CI     z  p-value
## Random effects model -0.3537 [-0.4882; -0.2192] -5.15 < 0.0001
## 
## Quantifying heterogeneity:
## tau^2 = 0.3666; H = 747523.89 [747522.42; 747525.37]; I^2 = 100.0% [100.0%; 100.0%]
## 
## Test of heterogeneity:
##                  Q d.f. p-value
##  90524299681441.44  162       0
## 
## Results for subgroups (random effects model):
##                                    k                     95%-CI
## intervention = vegetation         54  0.0434 [-0.0088;  0.0956]
## intervention = grazing exclusion  12  0.3481 [ 0.3412;  0.3550]
## intervention = soil               97 -0.4835 [-0.5740; -0.3929]
##                                                  Q   tau^2    I^2
## intervention = vegetation          775623742262.70  0.0075 100.0%
## intervention = grazing exclusion           6686.32 <0.0001  99.8%
## intervention = soil              11085228014467.45  0.1148 100.0%
## 
## Test for subgroup differences (random effects model):
##                       Q d.f.  p-value
## Between groups   446.83    2 < 0.0001
## 
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
#metabias(m)

#use random effects mean and var estimate to plot

#do t-tests here too

#outcomes
m <- metagen(lrr, var.es, studlab = ID, byvar = outcome, subset = paradigm == "active", comb.fixed=FALSE, data = mdata.all)
summary(m)
## Number of studies combined: k = 333
## 
##                                       95%-CI    z  p-value
## Random effects model 0.1947 [0.1033; 0.2862] 4.17 < 0.0001
## 
## Quantifying heterogeneity:
## tau^2 = 0.2649; H = 15484891.12 [15484889.86; 15484892.37]; I^2 = 100.0% [100.0%; 100.0%]
## 
## Test of heterogeneity:
##                     Q d.f. p-value
##  79607575146295024.00  332       0
## 
## Results for subgroups (random effects model):
##                     k                     95%-CI                    Q
## outcome = soil    143  0.1906 [ 0.0791;  0.3020] 79607530173242944.00
## outcome = plants  108  0.5161 [ 0.2783;  0.7538]          68146638.29
## outcome = animals  11  1.0512 [ 0.9693;  1.1331]                 0.05
## outcome = habitat  71 -0.4266 [-0.5776; -0.2755]        3395032696.46
##                    tau^2    I^2
## outcome = soil    0.2649 100.0%
## outcome = plants  0.3814 100.0%
## outcome = animals      0   0.0%
## outcome = habitat 0.0739 100.0%
## 
## Test for subgroup differences (random effects model):
##                       Q d.f.  p-value
## Between groups   344.41    3 < 0.0001
## 
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
#metabias(m)

m <- metagen(lrr, var.es, studlab = ID, byvar = outcome, subset = paradigm == "passive", comb.fixed=FALSE, data = mdata.all)
summary(m)
## Number of studies combined: k = 163
## 
##                                          95%-CI     z  p-value
## Random effects model -0.3537 [-0.4882; -0.2192] -5.15 < 0.0001
## 
## Quantifying heterogeneity:
## tau^2 = 0.3666; H = 747523.89 [747522.42; 747525.37]; I^2 = 100.0% [100.0%; 100.0%]
## 
## Test of heterogeneity:
##                  Q d.f. p-value
##  90524299681441.44  162       0
## 
## Results for subgroups (random effects model):
##                     k                     95%-CI                 Q   tau^2
## outcome = habitat  45  0.0091 [-0.0455;  0.0637]   775623742237.89  0.0075
## outcome = plants   21  0.3483 [ 0.3415;  0.3552]           6703.49 <0.0001
## outcome = soil     97 -0.4835 [-0.5740; -0.3929] 11085228014467.45  0.1148
##                      I^2
## outcome = habitat 100.0%
## outcome = plants   99.7%
## outcome = soil    100.0%
## 
## Test for subgroup differences (random effects model):
##                       Q d.f.  p-value
## Between groups   464.58    2 < 0.0001
## 
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
#metabias(m)

#t-tests for outcomes diff from 0 but you can see using the 95% CI what is different

#check metafor for interactions?? in these big models or are we ok??

#brainstorm on how to explore?? techniques


#save but cut all this.
#m.study <- metagen(lrr, var.es, studlab = study.ID, byvar = intervention, data = mdata)
#summary(m.study)
#funnel(m)
#metabias(m)
#forest(m, sortvar = intervention)

#aggregated data
#m1 <- metagen(lrr, var.es, studlab = intervention, data = simple.mdata)
#summary(m1)
#funnel(m1)
#metabias(m1)
#forest(m, sortvar = intervention)

#different var estimate
#m2 <- metagen(lrr, var.es, studlab = intervention, data = simple.mdata.var)
#summary(m2)
#funnel(m2)
#metabias(m2)
#forest(m, sortvar = intervention)


#m3 <- metagen(lrr, var.es, studlab = intervention, byvar = outcome, data = simple.mdata.2)
#summary(m3)
#funnel(m)
#radial(m3)
#metabias(m2)
#forest(m, sortvar = intervention)

#m4 <- metagen(lrr, var.es, studlab = intervention, byvar = outcome, data = simple.mdata2.var)
#summary(m4)
#funnel(m)
#radial(m4)
#metabias(m2)
#forest(m, sortvar = intervention)
#forest(m1, sortvar = intervention, fontsize = 10, xlim = c(-2, 2))

#forest(m2, sortvar = intervention, fontsize = 10, xlim = c(-2, 2))

#forest(m3, sortvar = intervention, fontsize = 10, xlim = c(-2, 2))

#forest(m4, sortvar = intervention, fontsize = 10, xlim = c(-2, 2))

#t-tests for lrr diff from 0
#mdata %>%
  #split(.$intervention) %>%
  #purrr::map(~sumz(.$lrr)) 
#effect sizes plots
#need better ci estimates
#ggplot(simple.mdata, aes(intervention, lrr)) +
 # ylim(c(-2,2)) +
#  geom_point(position = position_dodge(width = 0.5)) + 
 # labs(x = "", y = "lrr") +
 # coord_flip() +
 # geom_errorbar(aes(ymin=lrr-var.es, ymax=lrr+var.es), width=.05, position = position_dodge(width = 0.5)) +
 # geom_hline(yintercept = 0, colour="grey", linetype = "longdash")

#ggplot(simple.mdata.2, aes(intervention, lrr, color = outcome)) +
 # ylim(c(-2,2)) +
  #geom_point(position = position_dodge(width = 0.5)) + 
 # labs(x = "", y = "lrr", color = "") +
 # coord_flip() +
  #geom_errorbar(aes(ymin=lrr-var.es, ymax=lrr+var.es), width=.05, position = #position_dodge(width = 0.5)) +
  #geom_hline(yintercept = 0, colour="grey", linetype = "longdash")

#ggplot(mdata, aes(lrr, color = intervention)) +
  #geom_freqpoly(binwidth = 0.5, size = 2) + 
  #xlim(c(-10, 10)) +
  #labs(x = "lrr", y = "frequency", color = "") +
  #geom_vline(xintercept = 0, colour="grey", linetype = "longdash") +
  #scale_color_brewer()

#ggplot(mdata, aes(lrr, fill = intervention)) +
  #geom_dotplot(binwidth = 1) + 
 # xlim(c(-10, 10)) +
 # labs(x = "lrr", y = "frequency", fill = "") +
 # geom_vline(xintercept = 0, colour="grey", linetype = "longdash") +
 # scale_fill_brewer()

Interpretations

  1. Little to no replication of specific techniques.
  2. Intervention and outcome classification useful and appropriate.
  3. Active versus passive interesting and fascinating.
  4. Use only random effect model statistics because of heterogeneity and conceptual focus of data and synthesis.
  5. Active and passive are not the same - active is net positive and passive is net negative. Amazing finding.
  6. Then subset out each of active passive and run a meta. Found that both for active and passive, vegetation is a path forward for restoration (list techniques in paper) particularly active.
  7. Adding water does not generate significant effects, soil remediation is powerful active intervention and ignoring passive changes in soil is a dramatic and negative impediment to restoration.
  8. Habitat is challenging to actively restore but plants, animals, and soils are viable and significant positive outcomes that can be restored.
  9. Passive restoration confirms that habitat does not generally respond, plants can recover passively, but soil does not.
  10. Additional ideas - test ONE key moderator such as aridity index or length of experiment?